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I. INTRODUCTION 


It is generally recognized that the panel flutter phenomenon must be 
taken into account when designing the skin or exterior covering of high speed 
aircraft and missiles. Usable design criteria based on panel flutter aspects 
have not as yet been developed , and on the whole the designer is guided by 
previous experience and/or by experimental data obtained in special tests. It 
is clear, therefore, that unified design criteria are needed. 

Design criteria can be developed from quantitative data that suf- 
ficiently account for the effects of all pertinent parameters. At present, 
quantitative data for the panel flutter phenomenon are not available in such 
sufficient amounts. This lack of data is due in part to theoretical limita- 
tions, in part to the difficulty of performing good experiments, and in part 
simply to the fact that a large number of pertinent parameters are involved. 

A theoretical parametric study is therefore conducted here, having as its ob- 
jective quantitative flutter data for a moderately large range of the impor- 
tant parameters. It is hoped that these data will be useful both in the ulti- 
mate development of design criteria and in the intermediate task of evaluating 
theoretical techniques by means of suitable experiments. 

The theoretical parametric study is conducted for an array of simply 
supported, flat, rectangular panels, extending to infinity in the spanwise 
direction and to one panel in the chordwise direction. The low supersonic 
Mach number region is investigated, and the effects of Mach number, aspect 
ratio, panel material and altitude, are determined. The numerical results are 
presented in Figs. 2 to 29 * in the form of critical flutter boundaries in a 
generalized parameter plane. 


II, EQUATIONS OF MOTION 


Flutter equations of motion for an array of simply supported flat 
panels extending to infinity in the spanwise direction and to a finite number 
of panels in the chordwise direction are derived in [lj**. In the derivation, 
small deflection plate theory is used to represent the structural behavior, 
and since the low supersonic region is of interest, exact, linearized, three- 
dimensional aerodynamic theory is used to determine the aerodynamic pressure 
on the plate. An approximate solution to the resulting differential equations 


* Figures 1 through 37 appear in Appendix A. 

** Numbers in square bracket refer to the bibliography. 
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of motion for the transverse displacement w is obtained by the Galerkin 
method. The natural vibration modes of a continuous beam, simply supported at 
equally spaced intervals, are used to represent the chordwise and spanwise 
deflection shapes of the panel array. In the spanwise direction only the first 
mode is retained since it has been shown to yield the most critical flutter 
requirements. Thus, the transverse displacement of the panel array is approxi- 
mated by 


w (x,y,t) = Imag. 4 


N 


£ sin sin HI e* jt \ 


m=l 


where the generalized coordinates, q m , are complex in general. The flutter 
equations, finally obtained, are*: 


litM [(m 2 +s 2 ) 2 ^+(6 3 ffi 2 +6 2 s 2 ) M ]<^ - k 2 p,q~ 
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* Equations (l) follow directly from equations (5.7) of QlJ if the terms asso- 
ciated with large deflections are dropped and the coefficients regrouped. 
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and additional definitions are given in the list of symbols. 

The summation over the index r in the above terms results from an 
approximation introduced in the aerodynamic pressure term; namely, Bessel 
functions are approximated by a sum of circular functions (V) . The argument 
of the Bessel functions contains the chordwise coordinate x , and consequently 
the required number of terms, q , depends on the number of chordwise panels in 
the array. In the present computations for one chordwise panel, satisfactory 
accuracy is obtained using q = 6 . 


III. FLUTTER BOUNDARIES AMD METHOD OF COMPUTATION 


The flutter equations (l) are simultaneous algebraic equations for 
the generalized coordinates q^ • When N chordwise modes are used to repre- 
sent the deflection function, N equations are obtained, each equation being 
associated with a particular value of m • A nontrivial solution of the equa- 
tions marks a condition of sustained harmonic oscillation or neutral stability 
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for the panel array. Loci of such points plotted in a suitable parametric 
plane separate regions where flutter will occur from regions of no flutter, 
i.e., these loci may be interpreted as flutter boundaries. The numerical re- 
sults of the present study are presented in the form of such flutter bounda- 
ries . 


Before considering computation details it is useful to look briefly 
at the general form of the flutter equations. For convenience, the equations 
are first written in matrix form. 



( 2 ) 


where 


P - 
m,m 


.(L + M 

M 2 



A _ 

m,m 

B 


( 3 ) 


Nine basic parameters appear in the equations, namely, 

M - Mach number 

k - reduced frequency 

s - inverse of aspect ratio 

g - structural damping coefficient 

p - density in the free stream 

p 0 - material density 

c co - sound velocity in the free stream 

E - modulus of elasticity 

v - Poisson's ratio 

The Ajjj - are functions of M , k , and s . It is pertinent that the g 
are independent of the parameters associated with panel material and altitude, 
i.e., p , p , c w , E , and v . If and 62 are momentarily ignored 

(a x = Cy = 0) , it is seen that the material-altitude parameters appear only 
in a and p, . This means, in essence, that the nine basic parameters listed 
above, are consolidated into six fundamental groups, M , k , s , g , a , and 
p, . Thus, flutter boundaries can be obtained in general in terms of a and 
p, and it is not necessary to make separate computations for each combination 
of material and altitude of interest. This advantage can be extended to cases 
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of finite a x and by computing for an array of values 6^ and 6^ , 

and then interpreting the results in terms of O and o for various values 
of PgC,,^ . The procedure is somewhat indirect but the d?sadvantage is not 
serious . 


Returning no w to the details of computation it is recalled that 
flutter boundaries are to be obtained, consisting of loci of points represent- 
ing nontrivial solutions of the flutter equations (2). For (2) to have a 
nontrivial solution it is necessary that 


Det. 



(4) 


which is commonly called the flutter condition. The problem, then, is to de- 
termine those combinations of the parameters which satisfy the flutter con- 
dition. 


Since the determinant of (4) is complex and only real values of the 
parameters have physical significance, two parameters must be kept free in 
order to satisfy (4). It is convenient to take a and p, as the free param- 
eters. Flutter boundaries in the a-^ plane are then obtained for a selected 
set of values M , s , g , 6 1 , $ 2 by determining the a-ji pairs which sat- 
isfy (4) for various judiciously chosen values of reduced frequency k . Two 
procedures were used for finding these a-^ pairs. The first procedure to be 
described was used initially and then discarded in favor of the second for 
reasons that will be pointed out. 

The first method consists of a trial and error procedure for deter- 
mining the <x— (j, pairs which cause both the real and imaginary parts of the 
determinant to vanish simultaneously. For a particular value of a and a 
number of values of ^ the real part, R , and the imaginary part, I, of the 
flutter determinant are calculated. The values of y, for which R = 0 (R q ) 
and I = 0 (l Q ) are obtained by interpolation. Next, a is varied and the 
process repeated. In this manner the functions R Q (a) and I Q (a) iecome 
known. The intersections of R Q (a) and I Q (a) determine the flutter points 
and are obtained by interpolation. Although this method yields satisfactory 
results, it leaves something to be desired in the way of economy and con- 
venience of application. The principal difficulty arises from the fact that 
the functions R (o) and I (a) are multi-valued and conposed of several- 
branches in the a-p. plane. 
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The second method for determining and y, is as follows. The 
determinant of (4) can be written as a polynomial in (j, with complex coef- 
ficients that depend on a > 


- a oM- m + a lU m ' 1 + + a m 


Because D should vanish and p, must be real, the flutter condition becomes 
f(ji) - g, ra + (ap/ao^® -1 + + ( a m / a 0 )r = 0 


g(ki) = ( a ]/ a o V*"' 1 + + (^/sq)! * 0 


where the R and I subscripts mean "real part of" and "imaginary part of," 
respectively. Thus, at flutter, f(p,) and g(^) must have at least one common 
root. The condition that two polynomials in one variable have a common root 
can be obtained from Euclid's algorithm (*30, as follows. Since f is of 
higher degree in p, , it can be divided by g to give 


f(n) = Q 0 (u)g(n) + RjU) 

where Q q is the quotient polynomial and R^ is the remainder polynomial. 
Continuing the process by dividing g by R^ , etc., the following system of 
identities results. 

g(n) = Qp(n)Rp(|i) |jl ) 

R l(n) = Q 2 (h)R 2 (u) + R 3 (|i) 


^p-l^u) - Qp(g,)Rp(n) + Rp +1 
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Now g , R-j^ , R 2 


are polynomials of decreasing degrees in y, so that 


after a finite number of steps a remainder, Rp+p > I s reached which is in- 
dependent of y, . It can be shown that the condition that f and g have 


at least one common root is R ,, = 0. 

V . p + l 

pj = 0. Since a Q , a^ , 

function of o, and the flutter condition is R .(a) = 0. 
corresponding to flutter can thus be det ermined Dy the usual interpolation 
procedures. The associated values of (j, follow from Rp(u) = 0* 


Further, the common roots can be ob- 
, a^ depends on a , ^p+p * s 3 
The values of a 


This second method for determining a and y, is more straightfor- 
ward than the first because only one parameter, a , is varied and y follows 
as an interim result of the procedure. Also pertinent is the fact that the 
function Hp+p(a) , which must be monitored in the computer program, is single- 
valued in contrast with the multi-valued functions R Q (a) and I Q (a) of the 
first method. Most of the numerical results were obtained using this second 
method. 


Computation of the flutter boundaries are carried out on the IEM 
7090 computer. In all cases the chordwise deflection shape of the panel array 
is approximated by the first four natural vibration modes of a simply supported 
beam. As noted earlier, the spanwise deflection is approximated by the first 
mode of a continuous beam, simply supported at equally spaced intervals. The 
scope of the parametric study is discussed In the following section. 


IV. SCOPE AND RESULTS OF PARAMETRIC STUDY 
' - ' ■ ' ■— "1 ' ■ " 1 — 1 — 1 


Flutter characteristics are computed for semi-infinite panel arrays 
of the type shown in Fig. 1. The desire of the parametric study is to examine 
the effects of Mach number, length/width ratio, panel material, altitude, and 
initial membrane forces. Computations are therefore carried out for length/ 
width ratios, s , of l/4, l/2 , 1 , 2 , and 4 at each of six Mach numbers, 

1.1 , 1.25 , 1.35 , ^1T* , 1.5 , and 2,0.* Generalized parameters are em- 
ployed so that the results are applicable to a wide variety of panel materials 
and altitudes. At the least, the altitude range sea-level to 50,000 ft. is 
covered for brass, steel, titanium, aluminum, glass, magnesium, and copper. 
(The material-altitude coverage is described more precisely later. ) The 
effects of membrane forces are not examined because of time limitations; all 

results are obtained for zero membrane forces (a = a = 0). 

x y 


* Two cases, s = 4 and M = 1.1 and 1.25 were not completed because of 
numerical difficulties. 
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The numerical results of the parametric study are presented in Figs. 

2 through 29 and Tables I through V.* Before examining these results in detail, 
it is useful to first consider the parameters used for displaying the flutter 
boundaries. It may be recalled that, earlier, two general parameter groupings 
a and p, were identified, which consolidate the specific quantities associ- 
ated with panel material and altitude. Use is made here of a, and p to 
facilitate concise presentation of the numerical results. In the figures, the 

parameter on the ordinate is a 1 / 3 = t (E/pc^l-v 2 .] 1 / 3 . On the abscissa the 
parameter is (p'Va)"''/ 2 = jc^Pg/p] [jp s ( 1- which form was selected to 

eliminate the thickness ratio t • The two parameters are made up of the 
thickness ratio t ; material constants p g , E , and v , and airstream proper- 
ties p and c<n • Therefore, the computed stability boundaries can be applied 
conveniently, either to particular material-altitude combinations or to com- 
binations of material and wind tunnel conditions p and c* . The extent of 
coverage provided by the parametric study is indicated in Fig. 36 where 

(p^/a.) 1 / 2 is plotted against altitude for some typical panel materials. 

The ordinate parameter t [E/pc C0 2 (l-u 2 )J 1 / 3 is closely related to a 
dynamic pressure parameter, rfflE/q] 1 / 3 , which is commonly used in the litera- 
ture for correlation of experimental data. The two parameters are related by, 

tCe/pc^I-w 2 )] 1 / 3 x ^(l-v 2 )^ 2 ] 1 / 3 = t [jjE/q] 1 / 3 

Returning now to consideration of the numerical results, it is re- 
iterated that the results are completely contained in Figs. 2 through 29, 
where the computed flutter boundaries are shown in detail, one boundary to 
each figure.** In subsequent figures the boundaries are reproduced together, 
in lesser detail, for purposes of comparison. The region above the boundary 
is the stable region as indicated. The stability boundaries, therefore, de- 
fine the minimum thickness that will prevent flutter. A somewhat arbitrary 
terminology Is used to indicate the form of the flutter vector. A "first mode 
dominant" boundary, for instance, is one for which the flutter vector consists 
predominantly of the first generalized coordinate, q^^ , the magnitude of the 
second largest being less than q-jj/2 . A "first and second mode coupled" 
boundary, on the other hand, is one where q is largest and q second 
largest, with q 2 > q-j/2 . 1 2 


* Tables I through V are shown in Appendix B. 

** The numerical data from which the boundaries are constructed are given in 
Tables I to V. In the tables, the generalized coordinates, q m , are 
given in normalized polar form, referenced to the magnitude and phase of 
the largest coordinate. 
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Two interesting features are seen in the detailed stability bound- 
aries (Figs. 2 to 29). First of all, there are cases where the critical flutter 
vector changes abruptly from one mode of dominance to another. For instance, 
for M = ffsT and s = l/4 (Fig. 5), an aluminum panel will flutter most 
readily in the third mode at sea level, the fourth mode at higher altitudes, 
and finally in the first and second modes coupled at altitudes above approxi- 
mately 17,000 ft.' Brass, steel, or titanium panels on the other hand, will 
always flutter in the coupled mode. This switching of flutter modes, which is 
noted for all length-width ratios except s ■ 4 , is accompanied by significant 
changes in flutter frequency (see the tables). 

The second feature of interest is that the stability boundaries are 
only fairly flat in the high Mach number regime, i.e., the parameter 
TfE/pc^l-v 2 )] 1 / 3 at flutter is insensitive to changes in panel material and 
altitude only for M > Nf2~ , as theoretically expected. This point is clearly 
illustrated in Figs. 30 through 34. 

In Figs. 30 to 34 the stability boundaries for particular values of 
length-width ratio are superposed. The uppermost boundaries are critical for 
the Mach number range 1.1 to 2.0 considered in the parametric study. It is 
seen that the lower Mach numbers are critical for small s and that the criti- 
cal Mach number becomes larger with increasing s . The flutter vectors 
associated with the critical boundaries are as follows. 

s = l/4 ; first mode dominant 

s = l/2 ; second mode dominant 

s = 1 ; third mode dominant 

s = 2 ; fourth mode dominant 

s = 4 ; third and fourth mode coupled 

Thus, as s increases, the higher modes and higher Mach numbers become in- 
creasingly important. Also, examination of the tabulated data shows that, in 
general, coupling between modes increases as s increases. 

In Fig. 35 the uppermost boundaries from Figs. 30 to 34 are shown in 
superposition, to illustrate the pronounced stabilizing effect associated with 
increasing values of the length-width ratio s . The effect is greatest at 
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the smaller s values where the lower modes dominate, but the reduction in 
thickness required to prevent flutter is still appreciable in going from 
s = 2 to s = 4 . 


V. CONCLUSIONS AMD RECOMMENDATIONS 


The stability boundaries in Figs. 2 through 29 provide information 
on the flutter characteristics of flat panels under a wide variety of condi- 
tions. Detailed information for particular materials, altitudes, or wind 
tunnel conditions can be obtained conveniently using intermediate relationships 
of the type shown in Figs. 36 and 37. This detailed information should be of 
use to the designer in lieu of presently unavailable design criteria, and to 
the experimentalist in planning exp eriment s and analyzing exp eri mental observa- 
tions . 


A most critical need at this time is valid comparisons between 
theory and experiment . The data presented here should facilitate such com- 
parisons. 


For the future it is suggested first of all, that the length-width 
ratio s = 4 , treated here, should be further examined using six or more 
modes in the analysis. It is felt that the general characteristics cited in 
the previous section imply that the higher modes may have significant effect . 
It is further suggested that the effects of initial membrane forces and more 
than one panel in the chordwise direction be investigated, particularly since 
the theory is available (Eq.(l)). Finally, a parametric study similar to the 
one reported here, should be conducted for panels with partially clamped edges 
since it is expected that in practical application (and experiments) boundary 
conditions are somewhere between simply supported and clamped. 


Midwest Research Institute 
Kansas City, Missouri 
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APPENDIX A 


FIGURES 1 THROUGH 37 
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Fig. 2 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = 1.1, s * l/4, 
g = 0.01 (Table I) 
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Fig. 3 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = 1.25, s = l/4, 
g = 0.01 (Table I) 
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Fig. 4 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = 1.35, s = l/4, 
g = 0,01 (Table i) 
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Fig. 5 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = j s = l/^> 
g = 0.01 (Table i) 
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Fig. 6 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M - 1.5, s = l/4, 
g = 0.01 (Table I) 
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Fig. 7 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = 2.0, s = l/4, 
g = 0.01 (Table I) 
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Fig. 8 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = 1.1, s = l/2, 
g = 0.01 (Table II) 



Fig. 9 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = 1.25, s = l/2, 
g = 0.01 (Table II) 
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Fig. 10 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = 1.35, s = l/2, 
g * 0.01 (Table II ) 



Fig. 11 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = f*, s = 1/2, 
g = 0.01 (Table II) 
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Fig. 12 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = 1.5, s = l/2, 
g = 0.01 (Table II) 
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Fig. 13 - Critical Stability Boundary From a Four Mode Analyses 
for Simply Supported Flat Plates, M = 2.0, s = l/2, 
g = 0.01 (Table II ) 
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Fig. 14 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = 1.1, s * 1, 
g = 0.01 (Table III) 
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Fig. 15 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = 1.25, s = 1, 
g = 0.01 (Table III) 








Fig. 17 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = ^2", s = 1, 
g = 0.01 (Table III) 
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Fig. 18 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M « 1.5, s * 1, 
g * 0.01 (Table III) 
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Fig. 19 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = 2, s = 1, 
g = 0.01 ( Table III) 
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Fig. 20 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = 1.1, s = 2, 
g = 0.01 (Table IV ) 
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Fig. 21 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, Mj= 1.25, s = 2, 
g = 0.01 (Table IV) 
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Fig. 22 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M * 1.35, 
s * 2, g = 0.01 (Table IV) 
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Fig. 23 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M =«fi7 s = 2, 
g = 0.01 (Table IV) 
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Fig. 24 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = 1.5, s = 2, 
g = 0.01 (Table IV) 
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Fig. 25 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = 2, s = 2, 
g = 0.01 (Table TV ) 
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Fig. 26 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = 1.35, s = 4, 
g = 0.01 (Table V) 
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Fig. 27 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = f ^2~, s = 4, 
g = 0.01 (Table V) 
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Fig. 28 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = 1.5, s = 4, 
g = 0.01 (Table V) 
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Fig. 29 - Critical Stability Boundary From a Four Mode Analysis 
for Simply Supported Flat Plates, M = 2, s - 4, 
g = 0.01 (Table V) 
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Pig. 32 - Variation of Critical Stability Boundaries With 
Mach Humber for s = 1 
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Fig. 33 - Variation of Critical Stability Boundaries With 
Mach Number for s = 2 
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ALTITUDE -FT x 10' 



Fig. 36 - Variation of Parameter [ c coP s /pJ [p s (1-v 2 )/e] 1 / 2 
Tilth Material and Altitude 



Fig, 37 - Variation of Density Ratio p/p g With Panel 
Material and Altitude 
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TABLE I (Continued) 
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FLUTTER FREQUENCY AMD VECTOR FROM FOUR MODE ANALYSIS 
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